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Grassmann integrals by machine 

Michael Creutz a * 

a Physics Department, Brookhaven National Laboratory, P0 Box 5000, Upton, NY 11973-5000, USA 
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I present a numerical algorithm for direct evaluation of multiple Grassmann integrals. The approach is exact 
and suffers no Fermion sign problems. Memory requirements grow exponentially with the interaction range and 
the transverse size of the system. Low dimensional systems of order a thousand Grassmann variables can be 
evaluated on a workstation. 



In quantum field theory fermions are usually 
treated via integrals over anti-commuting Grass- 
mann variables [Q, providing an elegant frame- 
work for the formal establishment of Feynman 
perturbation theory. With non-perturbative ap- 
proaches, such as Monte Carlo studies on the lat- 
tice, these objects are more problematic. Es- 
sentially all approaches formally integrate the 
fcrmionic fields in terms of determinants depend- 
ing only on bosonic fields. When a background 
fermion density is present, as for baryon rich re- 
gions of heavy ion scattering, these determinants 
are not positive, making Monte Carlo evaluations 
tedious on any but the smallest systems [[|. This 
problem also appears in studies of many electron 
systems doped away from half filling. 

Here I explore the possibility of directly eval- 
uating the fermionic integrals, doing the neces- 
sary combinatorics on a computer [Q. This is in- 
evitably a rather tedious task, with the expected 
effort growing exponentially with volume. Nev- 
ertheless, in the presence of the sign problem, all 
other known algorithms are also exponential. My 
main result is that this growth can be controlled 
to a transverse section of the system. I illustrate 
the technique with a low dimensional system in- 
volving of order a thousand Grassmann variables. 

I begin with a set of anti-commuting Grass- 
mann variables {tp}, satisfying tpiipj + tpjtpi — 
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0. Integration is uniquely determined up to 
an overall normalization by requiring linearity 
and "translation" invariance J dtpf(tp + ip') = 
J dtpf(tp). I normalize things so that 



dtp ip = I 



dip 1 = 0. 



(1) 



Consider an arbitrary action S(tp) inserted into a 
path integral. I want to evaluate 



Z = dtp n . . . dtpi 



(2) 



Formally this requires expanding the exponent 
and keeping all terms containing exactly one fac- 
tor of each tpi. 

I first reduce the required expansion into op- 
erator manipulations in a Fock space. Intro- 
duce a fermionic creation-annihilation pair for 
each fermionic field, tpi <-> { 
isfy the usual relations [oj , at 
is built by applying creation operators to the vac- 
uum, which satisfies Oj|0) = 0. It is convenient 
to introduce the completely occupied "full" state 
\F) Eat ^ 

integral as the matrix element 



,t 



ai} 

= S. 



These sat- 
y . The space 



a ; ...a\a\ |0). Then I rewrite my basic path 



Z 



<0| 



\F)- (3) 
a non-vanishing contribution re- 



Expanding e 
quires one factor of aj for each Fermion. This is 
the same rule as for Grassmann integration. 

I now manipulate this expression towards a se- 
quential evaluation. Select a single variable tpi 
and define Si(a) as all terms from the action in- 
volving a factor of Oj. I define the complement Si 



2 



as anything else, so that S = Si + Si. I assume 
a bosonic action so that Si and Si commute and 
Z = (0\e Si e Si \F) . Since Si contains no factors 
of Oi, the occupation number for that variable, 
rii = a\a,i, vanishes between the two factors. I 
thus can insert a projection operator 1 — m 



Z = (0|e Si (1-ni) e Sl |.F} 



(4) 



Since 1 — rij projects out an empty state at loca- 
tion i, I trivially have a 2 ; (1— n,) = 0. I can replace 
Si with the full action. Also, since af = 0, the 
right hand factor expands as e Si = 1 + Si, giving 



Z= (0|e s (1-ni) (1 + 501^). 
Repetition gives my main result 



Z = ( 



Y[((i-m)(i + Si)) 



F 



(5) 



(G) 



assuming constants in the action are removed. 

Eq. (^|) summarizes the basic procedure. Cre- 
ate an associative array (hash table) to store gen- 
eral states of the Fock space. For a given state 
\fa) = J2 s Xs\s) store the numbers \s labeled by 
s) . Initially this table only contains the one entry 
for the full state. The algorithm then loops over 
the Grassmann variables. For a given fa, first ap- 
ply (1 + Si) to the stored state. Then empty the 
location with the projector 1 — rij. After all sites 
are integrated over, only the empty state survives, 
with the desired integral as its coefficient. 

The advantages appear with a local interaction. 
All sites previously visited are empty, and involve 
no information. Unvisited locations outside the 
interaction range are still filled, and also involve 
no storage. All relevant states are nontrivial only 
for unvisited sites within range of previously vis- 
ited sites. Sweeping through the system in a di- 
rection referred to as "longitudinal," we only need 
keep track of a "transverse" slice of the model. 
This is illustrated in Fig. ([!]). Although the di- 
mension of the Fock space is two to the number 
of Grassmann variables, the storage requirements 
only grow as two to the transverse volume. 

The approach is exact, with no sign problems. 
The complexity grows severely with interaction 
range, probably limiting practical applications to 
short range interactions in low dimensions. Note 
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Figure 1. Integrating out sequentially, visited 
sites are empty and out of range sites are filled. 
When integrating the site labeled fa only those 
sites labeled "X" are undetermined. 



that the effort only grows linearly with the lon- 
gitudinal dimension, allowing very long systems. 
This discussion has been in the context of "real" 
Grassmann variables. For "complex" variables 
treat ip and fa independently. 

In the transverse direction the boundary con- 
ditions are arbitrary, but longitudinal boundaries 
should not be periodic. To make them so re- 
quires maintaining information on both the top 
and bottom layers of the growing integration re- 
gion, squaring the difficulty. Note that the tech- 
nique is similar to the finite lattice method used 
for series expansions and closely related to a 
direct enumeration of fermionic world lines [^) . 

As a test, consider a spin-less fermion hopping 
along a line of sites. I introduce a complex Grass- 
mann variable on each site of a two dimensional 
lattice and study 



(dipdfa)e 



S t +S h +S! 



(7) 



with the various terms 

S t = Em V^Om -fa,t-i) 

Sh = k J2 lt fa t fa+i,t + fa+i, t fa,t (8) 

s i = sE llf ^tt^i,tfa+i,t^+i,t- 

I take Nt sites in the time t direction and Ni 
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F=log(Z)/V 

50 site chain, k=1 
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Figure 2. The free energy F = \og(Z)/NiNt with 
a four fermion interaction as described in the text, 
plotted as a function of the number of time slices 
N t . The chain has Ni = 50 sites. Points are 
shown for k — 1 and g = ±1. 



Figure 3. The occupancy of a 20 site chain as a 
function of the chemical potential scaled by the 
number of time slices. The filling occurs earlier 
or later depending on the sign of the coupling. 



spatial sites. The one-sided form of the tem- 
poral hopping insures an Hcrmitean temporal 
transfer matrix[||. This model Bosonizes into an 
anisotropic quantum Heisenberg model, a fact not 
being used here. 

I treat t as my "transverse" coordinate, growing 
the lattice along the spatial chain. Fig. (||) shows 
the Nt dependence for the free energy with k = 1 
and g = ±1. Here memory requirements were 
reduced by using time translation invariance after 
integrating each layer. The Nt = 14 points were 
run on the FJKEN/BNL Supercomputer. 

With Monte Carlo methods, a chemical poten- 
tial term can be highly problematic due to can- 
celations. Here, however, it is just another local 
interaction of negligible cost. As an illustration, 
take S — St + Sh + Si + Sm with 



s M = Mj2i>:, t ^. 



(9) 



i.t 



This regulates the "filling," which can be approx- 
imately monitored as (1 + M)-r£. I include the 
factor of 1 + M to compensate partially for fi- 
nite Nt artifacts. Fig. (|J) shows the filling as a 
function of MN t on an N t = 8 by Ni = 20 lat- 
tice with a spatial hopping parameter of k = 0.1. 
Here I made a crude extrapolation in chain length 
by defining F = log ^ z ^n^-i) ) ■ Note how the 



four fermion coupling enhances the filling. 

An obvious system for future study is the Hub- 
bard modclQ. This requires 4 Grassmann vari- 
ables per site corresponding to tp* and ip for spins 
up and down. Higher spatial dimensions strongly 
increase the size of the transverse volume and will 
limit practical system volumes, but this may be 
compensated for by the lack of sign problems. 
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